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Abstract 

Amorphous solids can flow over very long periods of time. Solid flow can also be artificially 

enhanced by creating defects, as by Ion Beam Sputtering (IBS) in which collimated ions with 

energies in the 0.1 to 10 keV range impact a solid target, eroding its surface and inducing formation 

of nanometric structures. Recent experiments have challenged knowledge accumulated during the 

last two decades so that a basic understanding of self-organized nano-pattern formation under IBS 
O ' 

is still lacking. We show that considering the irradiated solid to flow like a highly viscous liquids 

can account for the complex IBS morphological phase diagram, relegating erosion to a subsidiary 

>. 

role and demonstrating a controllable instance of solid flow at the nanoscale. This new perspective 



o 
o 



X 



can allow for a full harnessing of this bottom-up route to nanostructuring. 



The difference between solids and liquids is straightforward for the pedestrian: fluids flow 
and solids don't. This distinction is not so clear cut when tryingto be more quantitative. 



a 



For instance, when observed over long periods of time, glaciers [l| or lead pipes also flow 
in spite of being hard solids. Solids can indeed flow when observed at the proper time 
scale, the main cause being defect rearrangement. Unfortunately, experimental observation 
is difficult for macroscopic systems because of the long time periods involved (decades or 
even centuries). 

Solid flow at the nanoscale can be also observed or even enhanced at short times by 
creating artificially a large number of defects (see an sketch in Fig. [T]). Using ions with 
average kinetic energy E from, say, 100 eV up to 10 keV, this can be achieved by Ion Beam 
Sputtering (IBS) of solids that become amorphous under this type of irradiation (such as 
semiconductors) [2]. The ions create permanent defects in a surface layer whose thickness 
is of the order of the average penetration depth, trigger material ejection (sputtering) from 
its surface and induce material redistribution and flow 

Among the intriguing features of IBS, its efficiency to produce surface nanopatterns over 
large areas on top of a wide variety of targets, with a high degree of reproducibility p, |6j, 
furnishes the technique with a large potential for applications. Notably, samples that were 
unaffected by contamination issues have not been available until recently 3, |8[, a fact that 
has hindered a correct understanding of the mechanisms controlling this nanostructuring 
procedure. The ensuing experiments on Si have challenged the current paradigm according 
to which pattern formation is controlled by the competition between erosion processes and 
atom diffusion on top of the target (the so-called Bradley-Harper (BH) mechanism J9|j. In 



fact, none of the current phenomeno logical models that elaborate on this view 10Ml2j] are 
able to capture the complex morphological features elucidated by the experiments of Madi 
and coworkers, including the lack of pattern formation at small incidence angle 9 or the 
existence of more than one type of morphological transitions in a (6, E) parameter space. 
In these theories, the target surface is described through an effective equation in which the 
various physical mechanisms are combined ad-hoc, defect dynamics being neglected. The 
same limitations occur even for alternative models in which the evolution of the surface 
height is explicitly coupled to that of the density of defects whose transport is confined to 
a thin surface layer, akin to pattern formation on e.g. aeolian sand dunes, see e.g. |13| and 
references therein. 



Here, we adopt a new approach to the problem by considering that the external amor- 
phous layer produced on the target has a finite thickness and can actually be described as 
a highly viscous fluid. Our theory relegates the role of erosion as the dominant pattern 
forming mechanism and describes how solid flow can lead to both smoothening and nanos- 
tructuring, and the morphological transitions among them. This fluid-like behavior of the 
solid thus accounts for the pattern forming properties recently observed, our theory being 
validated through comparison with available experimental data. In retrospect, it is inter- 
esting to note, as in the first experimental reports half a century ago |l4|, that many of the 
appealing patterns produced by IBS resemble waves on the surfaces of water or sandy dunes 



(see Fig. 3 in |15{|). These shapes are manifestations of the physical processes occurring in 
the material media beneath. Thus, IBS appears as an ideal context in which a connection 
between fluid and solid flow can indeed be established. 

In the energy range considered, an impinging ion interacts mostly with the nuclear struc- 
ture of the target, generating vacancies and interstitials. Both ions and defects can diffuse at 
different rates, annihilate by recombination or form clusters. As borne out from Molecular 
. Dynamics (MD) simulations, the overall effect of these processes is threefold: amorphizing 
16 1 and stressing the target [17|, and displacin g th e mean atom position inside the material 
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181 ]. Moreover, as also seen in experiments |2|, |l7|, ]l9|, the amorphous layer forms rapidly, 
its thickness and other mechanical properties becoming stationary after a fluence of order 
10 14 ions cm~ 2 (corresponding to a few seconds of irradiation for typical ion fluxes), the 
system evolves quasi-statically, and the density, p, of the target in the amorphous phase 
remains almost a constant. In contrast with the macroscopic solid flow of a glacier, ion 
fluxes employed in standard IBS experiments imply a slow driving that induces a time scale 
separation between atomistic relaxation rates (~ 1 ps _1 for e.g. collision cascades, adatom 
hopping attempts), and the ~ 10 s times in which significant morphological changes occur. 
This fact legitimates a continuum description. 

In the language of fluid mechanics, we assume amorphous layer flow to be incompressible, 
V • V = 0, where V = u\ + wh is the velocity field. Moreover, it can be considered as highly 
viscous because the radiation induced viscosity is around fifteen orders of magnitude larger 
than the viscosity of water [3|. Hence, Stokes flow applies, 

pDV/Dt = = b + V ■ T, (1) 



where T t j = —PSij + /i (diUj + djui) are the components of the stress tensor T, with P 
being pressure and /i viscosity. Eq. ([1]) is supplemented with proper boundary conditions 
involving the erosive terms and the stress at the interfaces. Specifically, a no-slip condition 
is assumed at the crystalline-amorphous interface and erosion and amorphization rates are 
taken equal as to guarantee stationarity of the amorphous layer thickness, see |l5j |. We factor 
the body force term, b, into a momentum exchange term per unit volume, times an angular 
contribution, b — |b| = $e^{9 — l)- The coefficient Je, which has units of a stress gradient, 
contains the coarse-grained information about the effect of the residual stress created in the 
target, due to ion-induced mass redistribution. Thus, Je can in principle be estimated from 
macroscopic measurements of stress (see below). Regarding the geometrical factor in b, it 
only depends on the difference between the angle of incidence 9 and the local surface slope 7 
(see Fig. [Q. Symmetry requires \l/ to be an even function of its argument. Additionally, in 
the limit of completely glancing incidence the ion does not collide with the target, so that \1/ 
must vanish when 9 — 7 = tt/2. We choose \1/ = cos(9 — 7) as the simplest function fulfilling 
these constraints. Moreover, a simple cosine has the meaning of a geometrical correction to 
the average ion flux. 

In order to find the conditions under which the system evolves into patterns or produces 
smooth surfaces, we consider an initially flat surface and determine the stability of such a 
solution to periodic modulations of spatial wavelength A = 2ir/q (q is the space wavenumber), 
whose amplitude evolves in time as e Uqt . All the information about stability is contained in 
the dispersion relation, u q . Hence, a modulation with wavelength A will grow provided the 
real part of the dispersion relation, cu' , is positive, and will be damped out for a negative u' q . 
The experimental surface will be flat only if uj' < for all wave vectors. The imaginary part 
of the dispersion relation, u'' contains information about how the pattern travels across the 

n 

surface. In our case, we find [15] 

(Je cos(29) + q 2 a) (-2hq + smh(2hq)) _, 
Ug 2qfx(l + 2h 2 q 2 + cosh(2hq)) q) U 

where fi is the ion induced viscosity and a is the surface energy. For a general angular 
function ^(9), the dispersion relation has the same form as in Eq. (T5]) but replacing cos 29 
with 80(^(9) SU16 1 ). The term uj' is the contribution to the dispersion relation coming from 
purely erosive terms. Eq. ([2]) has some interesting limits in Fluid Mechanics. For instance, 
for f E — > 0, it reduces to the classic result by Orchard 20( for the flow of a viscous layer of 



arbitrary thickness on top of an immobile substrate. The q — > oo limit of the latter reduces 
to Muffins' 21J cfassic rate of flattening for interface features when relaxation occurs through 



viscous flow u q = — crg/2/x, already invoked in the context of IBS experiments 22|, |23 ] . 



If we consider solid flow as the only pattern forming mechanism, we can drop the erosive 
contribution UJ' and re-evaluate this assumption a-posteriori after comparing with experi- 
ments. For the typical orders of magnitude of the layer thickness h (about a few nanometers), 
the sign of u' is determined by the sign of — cos(2#). Thus, Eq. (J2J) predicts a phase transi- 
tion at 9 C = 45° separating smooth surfaces (8 < 45°) from nanopattern formation (9 > 45°). 
The patterns generated in the unstable region have a characteristic wavelength given by 



A n = 27Tx 
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(3) 



-f E cos(26)- 

To be more quantitative, we need to determine the main parameters of the theory. As men- 
tioned above, the new parameter, /#, is related to the stress. There is a large experimental 
uncertainty in the measurement of stress, a property that depends moreover on the energy of 



the incident ions 
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25]. Different aut 



Si, in the range 200 MPa to 1.62 GPa 



lors report different values for ion-induced stress on 



171 ] . We are interested in the order of magnitude, so 



we take the geometrical mean (which accounts for the average order of magnitude) of both 
extreme values, thus, approximately 569 MPa. The body force is Je = stress/depth. Thus, 
the latter value of stress corresponds to /# = 0.424 kg nm~ 2 s -2 for a layer of 3 nm depth. 
In Fig. 2b, this value corresponds to the blue line, while the limiting values of 200 MPa and 
1.62 GPa are used to plot the red lines, which we take as our confidence interval. More ac- 
curate measurements of stress would provide improved estimations. Additional parameters 
have not been measured experimentally yet, as we are proposing a new approach. However, 
we can provide indirect estimates for their values when these are not explicitly available, as 
shown in Table 1 in 15] . With these parameters at hand, Eq. ([3]) and experimental data 
are seen to agree quite closely in Fig. [2b. 

Eq. ([3]) actually predicts power law behavior in the vicinity of the critical angle, A c ~ 
2-K{aj Je) 1 ^ 2 {Q — c )~ l l 2 nm, and u' q ~ f E h 3 (3fiay 1 (9 — 6 C ) 2 , which are the same functional 
forms of deviation from 8 C as experimentally reported [8[. In the language of Pattern For- 
mation [Tj, |26| , this morphological transition is a type II bifurcation, akin to an equilibrium 
second order phase transition, in which the typical wavelength diverges right at the tran- 
sition point {6 = 6 C ). Specifically, it is a continuous transition in the sense that the wave 



vector characterizing the pattern on the large angle side of the phase boundary emerges 
continuously from the origin for increasing 9, as seen on Fig. [2b. 

An additional prediction of our theory is that the critical angle does not depend on 
energy, viscosity, surface tension, or on the thickness of the amorphous layer (which may 
contain additional angle dependences). Actually, a clear transition between 45° and 50° has 
been reported experimentally by other authors 27|, |28| for other targets and/or irradiation 
conditions. The sharp theoretical transition corresponds to the purple dashed vertical line 
in Fig. |2^, emphasizing the generality of our predictions in terms of parameter dependences. 
Thus, beyond the particular case of Si targets, our theory provides a new striking result: both 
the morphological transition at the critical angle and the patterns observed for more oblique 
incidences are exclusively caused by solid flow and not by erosion-induced instabilities. This 
can be better understood through an analogy: the IBS system is similar to a liquid flowing 
down an inclined plane [29| with inclination angle 9 (see Fig. 2 in 15 1). The only difference 
between both systems is that, here, the external body force, b, is not constant (as gravity is in 
the fluid analog) but depends, rather, on the local inclination angle (through the dependence 
of $ on 7). Thus, in the IBS case the effective gravity is larger in those locations where 
the surface faces the external beam. Precisely, for 9 > 9 C , the net difference between the 
force on both sides of an undulation makes it shrink laterally and its amplitude grows due 
to incompressibility, this being the physical origin of the instability. 

Going back to the experiments summarized in Fig. [21 there are additional important 
morphological features in the (9, E) diagram besides the transition we have just discussed. 
Namely, there is another transition at small angles and low energies (indicated by the oblique 
green dashed line in the figure). This transition is of type I, namely, the typical wave leng th 
remains finite at the transition, akin to an equilibrium first order phase transition |7|, |26] in 
which the wave vector characterizing the system jumps disco ntinuously from q = to a finite 
value. Given that the solid flow mechanism is morphologically stable below 9 C , a destabilizing 
mechanism must prevail at low incidence angles in order to observe a patterned structure. A 
natural candidate is erosion, so far assumed to be the only cause of morphological instabilities 
in these systems [6|. We thus need to quantify the balance between the erosion and flow 
mechanisms by choosing a specific form of the erosive part of the dispersion relation, cJ' . Its 



classic form is 



6] 



K = -v E H {d)q 2 -JCq\ (4) 
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where z/f (9) is negative at small incidence angles. Given this choice, we show in Fig. [2H 
the shape of the full dispersion relation u' for increasing energies at a fixed small incidence 
angle, in such a way that the experimental phase boundary is crossed. Indeed, the less stable 
wave vector jumps discontinuously to zero from a finite value at the positive maximum of 
u' as expected for a type I transition. Also, the length scale of the pattern on the low 
energy side of the transition is in good agreement with the experimental one (Fig. [2fe), and, 
as the latter, shows a weak dependence on the ion energy. 

Note that the BH estimate that we are using for v§ H {6) is known from simulations, using 
e.g. the binary collision approximation, not to be very accurate numerically [6J. This may 
be the reason why the shape of the transition curve (purple triangles in Fig. [2^i) does not 
agree with the experimental one. However, the main conclusion is that BH type mechanisms 
can still compete with solid flow and induce instabilities at small incidence angles, provided 
the average ion energy is sufficiently small. 

Returning to the main questions considered in the introduction, ion induced solid flow 
allows to understand the intricate variety of morphologies observed in the laboratory at low 
to intermediate ion energies, and the transitions among them. In this way hydrodynamic 
interactions may constitute the previously hypothesized non-local effects % ll2|] required to 
account in particular for the observed type I transition. Interestingly, to some extent this 
fact parallels the case at higher ion energies (~ 1 MeV), although strong differences persist, 
like the predominance for the latter of electronic scattering over the nuclear contribution to 
the stopping power, leading to strong many body effects over binary collision cascades, and 



the ensuing viscoelastic flow [30] 



At any rate, one can exploit the analogy with fluid flow to translate both the theoreti- 
cal toolbox and the physical insight to understand more profoundly the complex interplay 
between relaxation mechanisms in IBS. Intuitively, the vision of the solid as a slowly mov- 
ing fluid on whose surface nanowaves can be created, can drastically impact the controlled 
production of nanostructures in semiconductors. Yet, our approach also raises many funda- 
mental questions and predictions whose assessment requires new atomistic simulations and 
experiments. Thus, many details need to be quantitatively determined, in a similar fashion 
as viscosity had to be tabulated in the past, following application of Newton's laws to Fluid 
Dynamics. This line of research opens new theoretical and practical issues to be elucidated 
in the coming years. 



Appendix A: Constitutive equations and boundary conditions 

The equations for the conservation of mass, 

V-V = 0, (Al) 

and momentum 

= b + V • T, (A2) 

introduced in the main text, need to be supplemented with proper constitutive equations 
and boundary conditions. We will assume the solid to flow as a Newtonian fluid, thus the 
stress tensor is given by 

T i: j = -PSij + fj, (diUj + djUi) , (A3) 

with i,j = {x,y}, where P is hydrostatic pressure (that is an output of the solution of 
the full problem) and /i is viscosity. We consider ions that are impigning in planes that 
are parallel to the (x, y) plane, and assume translational invariance in the perpendicular 
direction for simplicity. Given the fast stabilization of the thickness of the amorphous layer 



seen in experiments [2] and in MD simulations [171 . Il9| , we assume that the amorphous layer 
is in a stationary state in which it is bounded between two interfaces. Denoting by h^ 
the position of the free target surface above a reference plane and by h^ c ' the position of 
the crystalline-amorphous interface that exists beneath (so that h(x, t) = h^ — h^ in Fig. 
1 of the main text), we assume a no-slip boundary condition at the crystalline- amorphous 
interface, while we request a balance of stress at the free interface along the local normal 
direction 

l + (h x a) y " (l + (/^)2)3/2 ' 

and in the tangential direction 

{Tzz ~ T xx )hj a) + (1 - {h^f)T xz da | t 
1 + (hi a) y ~ ds 

Here, II and r are external normal and shear stresses applied directly on the surface, that 
will be henceforth neglected, a is the surface tension of the amorphous-vacuum interface, 
and s is arc-length. Explicit non-zero values for these stresses could be used for further 
improvements of the theory — for instance, one could consider contributions to II and r that 
are induced by implanted ions in the near-surface region. 
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Finally, the kinematic condition that the material derivative is zero at each interface leads 
to evolution equations for both h^ and h^ as 

w-d t h^-ud x h^= Jam , (A6) 

w - dth® - ud x hP> = j er , (A7) 

where j er and j am are rates of erosion and amorphization, respectively. For instance, the 
erosive (BH-type) dispersion relation enters directly through j er , so that its contribution 
to the full dispersion relation is additive, as used in the main text. Again, due to the 
assumed stationarity of the amorphous layer thickness, we will take both rates to be equal 
on average. Moreover, for simplicity, we assume that the crystalline-amorphous interface is 
flat, since it can be shown that deviations from flatness only influence the imaginary part of 
the dispersion relation, while only the real part influences ripple formation and its stability. 
Finally, the linear dispersion relation is obtained by solving the linearized version of the 
previous set of equations for a periodic perturbation of the flat solution, h(x,t) = e Uqt ~ iqx , 



where h = h,( a > — hS c ' following the standard procedure (for details, see e.g. Ref. 26|). 



Appendix B: Predictions of the theory in the proximities of the type II transition 

The typical length scale of the pattern above 9 C = 45° is given by 



Expanding the cosine in the proximity of the critical angle, we find 

cos(2fl) = cos(2fl c ) - 2 sin(26 c )(6 - C ) + • • • ~ -2(0 - 9 C ), (B2) 

so that 

"~ } -0c)- 1/2 , (B3) 

as observed in the experiments of Ref. [8|. Using this approximation in the full dispersion 
relation and assuming that h/X Cl we find that 

4 - f^ - C )\ (B4) 

also in agreement with the experiments. Eqs. (10)-(11) are as quoted in the main text. As 
stressed there, the critical angle depends only on the form of *£?(9 — 7) so that the numerical 



values of h, /i, a, etc. do not affect the predicted value of 9 C . In Fig. [3]we show, schematically, 
how the ensuing surface morphology can be tuned by changing the incidence angle. 

Appendix C: Numerical values of the parameters 



Parameter 


Value 


Reference 


Surface tension (a) 


1.34 J/m 2 


[34] 


Amorphous layer thickness (250 eV) 


3 nm 


&M 


Viscosity (/u) 


6 x 10 8 Pa s 


[35]* 


Stress range 


200 MPa-1.62 GPa 


[8,12 



.a 



TABLE I. Numerical values of the parameters used in the main text. *Mayr et al. |3J report the 
value \i = 10 9 Pa s, which also leads to good agreement with the results presented in this work. 



The type I transition at low angles and small energies arises from the competition between 
the stabilizing solid flow (through parameter /#) and a destabilizing mechanism. Until 
more refined measurements can be done at these incidence angles, we postulate erosion 
(sputtering) as described by BH as a natural candidate. We obtain estimates of the erosion 
rates from the TRIM software [31|. Thus, for Si at normal incidence 32[, ve = 3.6 nm 2 s _1 . 
The dependence of this value with energy is linear [10| . Similarly, the amorphous layer width 
scales with energy as h ~ E 2m with m ~ 0.31 [28J for the range of energies considered. On 
the other hand, /# is related to stress and to the amorphous layer thickness. In Ref. [24] 
the connection between stress, E, and energy is given by 

£1/2 



S oc 



r + kE 5 / 3 ' 



(CI) 



where r = 1.7 and k = 0.0002 are taken from [33[ for C (these parameters seem not to be 
available for Si). With these values, we use the constant /C in Eq. (4) of the main text as a 
free parameter, finding /C = 38.2 nm 4 s _1 in order to have a characteristic wavelength around 
75 nm. Estimates for other parameters employed in the main text appear on Table 1. 

The lack of precise parameter values for Si introduces uncertainties in the quantitative 
value of the dispersion relation. However, the qualitative form of the full linear dispersion 
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relation u' seems to be valid, and stresses the non-local character of our theory as a main in- 
gredient to account for the Type I transition. We expect further developments in simulation 
techniques as well as more accurate experiments to improve the present predictions. 

Appendix D: Analogy with a fluid flowing down an inclined plane 

In Fig. H] we provide a cartoon to illustrate the analogy between ion induced solid flow 
and a genuine fluid flowing down an inclined plane. The two systems are mathematically 
equivalent once a rotation of angle 9 is performed in the laboratory reference frame. The 
main difference between both systems comes from the corresponding gravity fields. While 
for the standard fluid problem gravity is a constant, for the erosion system the effective 
gravity changes from point to point on the surface, its value depending strongly on the 
local surface orientation with respect to the ion beam. Recalling the assumption on fluid 
incompressibility that has been made, this fact accounts physically for the morphological 
instability that appears at 45°, as discussed in the main text. 
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FIG. 1. Schematic view of an IBS experiment. An incident energetic ion (red) impinges onto the 
target at an angle 9, inducing a collision cascade that amorphizes a region of width h(x,t) through 
the creation of vacancies and interstitials. At long time scales the amorphous solid flows as a highly 
viscous fluid. Here 7 is the local slope of the surface profile. 
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FIG. 2. a) Experimental observations for Ar + irradiation of Si targets, taken from [7|, |8|. Red 
and blue: ripple formation; black: flat surfaces; purple circles: ripple formation; purple crosses: 
hole formation. The dashed green lines indicate the experimental phase boundaries for a type II 
phase transition (vertical line at = 48°) and for a type I transition (oblique line). The present 
theoretical predictions for these boundaries are shown as the purple dashed line at 6 = 9^, heo = 45°, 
and as the purple triangles, respectively, b) A c (0) at the type II transition. Experimental data 



from 



7, 



8|] are shown as black symbols, the green line is the experimental phase boundary, and the 
blue line is the present theoretical prediction, bounded within the red lines (confidence interval), 
c) Linear dispersion relation at type II transition, d) Linear dispersion relation at type I transition 
(note how the maximum of the curve gets positive at a finite value of q). e) Different typical 
wavelengths are obtained for different energies on the low energy side of this transition: Circles, 



experimental data from 



crosses: present theoretical prediction. 
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FIG. 3. The surface morphology changes with the incidence angle. This cartoon illustrates the 
transition at t ,! ieor = 45° from flat to patterned surfaces. The solid black arrow indicates the 
direction of the ion beam. Above the critical angle, the characteristic size of the pattern is smaller 
for larger angles, according to A c = 2ir[2a/(— /e COS20)] 1 ' 2 . 
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FIG. 4. Analogy between IBS and fluid flow down a inclined plane. By rotating the frame of 
reference so that the ion impacts vertically, the IBS system can be seen as fluid flow down an 
inclined plane. The body force b is equivalent to a non-homogeneous gravity field that depends 
on the local slope value at the free interface. 
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FIG. 5. Water ripples (top) of characteristic wavelength 1 cm vs Ion Beam Sputtering (IBS) ripples 
(bottom) on Si. Water ripples were obtained by shaking a photography bucket filled with water 
(courtesy of R. Vida, Universidad Pontificia Comillas). Silicon ripples were obtained by irradiation 
with Ar + . The width of the observation window is 1 /im 2 (courtesy of L. Vazquez, Instituto de 
Ciencia de Materiales de Madrid-CSIC). The naked eye similarity between these macroscopic and 



microscopic patterns were already noted in early works on IBS 



uh- 
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